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Testing the Gaussianity of the COBE-DMR data with 
spherical wavelets 
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ABSTRACT 

We investigate the Gaussianity of the 4-year COBE-DMR data (in HEALPix pixeh- 
sation) using an analysis based on spherical Haar wavelets. We use all the pixels lying 
outside the Galactic cut and compute the skewness, kurtosis and scale-scale correla- 
tion spectra for the wavelet coefficients at each scale. We also take into account the 
sensitivity of the method to the orientation of the input signal. We find a detection 
of non-Gaussianity at > 99 per cent level in just one of our statistics. Taking into 
account the total number of statistics computed, we estimate that the probability of 
obtaining such a detection by chance for an underlying Gaussian field is 0.69. There- 
fore, we conclude that the spherical wavelet technique shows no strong evidence of 
non-Gaussianity in the COBE-DMR data. 

Key words: methods: data analysis-cosmic microwave background. 



1 INTRODUCTION 



The cosmic microwave background (CMB) provides a unique 
tool for investigating the formation of structure in the Uni- 
verse. In particular, studying the Gaussianity of the CMB 
temperature fluctuations allows us to distinguish between 
two competing theories of structure formation: the standard 
inflationary model that predicts Gaussian fluctuations and 
topological defects that give rise to non-Gaussian signatures 
in the CMB. In order to test the Gaussianity of the CMB, a 
large number of methods have already been proposed in the 
literature (see e.g. Barreiro 2000 and references therein). In 
particular, many of them have been applied to the 4-year 
COBE-DMR data, although most have not yielded any de- 
tection of non-Gaussianity (see e.g. Mukherjee, Hobson & 
Lasenby 2000). 

Nevertheless, three apparently robust detections of non- 
Gaussianity in the 4-year COBE-DMR data have recently 
been reported. Ferreira, Magueijo & Gorski (1998) studied 
the distribution of an estimator for the normalised bispec- 
trum, finding that Gaussianity is ruled out at the 98 per 
cent confidence level with the non-Gaussian signal mainly 
concentrated on the multipole I — 16. Although this non- 
Gaussian signal is certainly present in the data, Banday, 
Zaroubi & Gorski (1999) have shown in a recent work that 
this non-Gaussian signal is most probably due to a system- 



atic artefact. However, Magueijo (1999) has found a new 
non-Gaussian signal above the 97 per cent confidence level 
using an extended bispectrum analysis, which is present even 
when the artefacts found by Banday et al. are removed. 

Complementary to the above bispectrum analysis, 
Pando, Valls-Gabaud & Fang (1998) applied a technique 
based on the discrete wavelet transform to Face and Face 
5 of the QuadCube pixelisation of COBE-DMR data. On 
computing the scale-scale correlations of the wavelet coef- 
ficients in certain domains of the wavelet transform, they 
found a significant non-Gaussian signal at the 99 per cent 
confidence level in Face corresponding to scales of 11 — 22 
degrees. However, no significant deviation from Gaussianity 
was found using the skewness and kurtosis of the wavelet 
coefficients at each of the considered domains. Mukherjee, 
Hobson & Lasenby (2000) (hereinafter MHL) have recently 
revised the previous work, taking into account that a large 
number of the computed statistics show no evidence of non- 
Gaussianity, and pointing out that the results depend crit- 
ically on the orientation of the data. They find that Gaus- 
sianity can only be ruled out at the 41 per cent confidence 
level in the DSMB data and at the 72 per cent level in the 
53-1-90 GHz coadded data and therefore that this analysis 
does not provide strong evidence for non-Gaussianity in the 
COBE-DMR data. 
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The above wavelet analyses were performed by apply- 
ing planar (Daubechies) wavelets to Face and Face 5 of the 
COBE-DMR QuadCube pixelisation in Galactic coordinates 
(i.e. the faces centred on the North and South Galactic poles 
respectively). This procedure thus uses only one-third of the 
COBE-DMR data and, furthermore, can lead to distortions 
of the CMB fluctuations when moving from the sphere to 
the planar faces of the QuadCube. Therefore, in this pa- 
per, we use orthogonal spherical Haar wavelets (SHW) (Gi- 
rardi & Sweldens 1995, Sweldens 1995), which are better 
suited to analysing data over large regions of the sky (such 
as the COBE-DMR maps). These wavelets were introduced 
as a generalisation of planar Haar wavelets to more gen- 
eral spaces than R" . Indeed, Tenorio et al. (1999) have used 
SHW as a tool for studying the spatial structure, denoising 
and compression of CMB maps. We use SHW to perform 
a similar analysis to those of Pando et al and MHL, based 
on the skewness, kurtosis and scale-scale correlation of the 
wavelet coefficients, in order to search for evidence of non- 
Gaussianity in the 4-year COBE-DMR data. A hierarchical 
pixelisation scheme which is particularly well-suited to the 
application of such a wavelet decomposition is HEALPix 
(Gorski, Hivon & Wandelt 1999). Therefore, we apply the 
former analysis to the COBE-DMR data in HEALPix pixeli- 
sation (Banday et al. 2000). We also include a Galactic cut 
derived by a simple propagation of the customised Galac- 
tic cut of Banday et al. (1997) and Bennett et al. (1996) 
to the HEALPix pixelisation. Owing to the characteristics 
of the SHW transform, it is straightforward to use all the 
pixels lying outside the Galactic cut, which constitutes ap- 
proximately two-thirds of the total number of COBE-DMR 
pixels. 



2 THE WAVELET ANALYSIS 
2.1 Spherical Haar Wavelets 

Orthogonal SHW are an example of the so-called 'second 
generation wavelets' (Sweldens 1995, Schroeder & Sweldens 
1995). These wavelets are not dilations and translations of a 
given function and can therefore be adapted to more general 
spaces than but, at the same time, they still enjoy all 
the useful properties of planar wavelets, such as good space- 
frequency localisation and a fast transform algorithm. 

The temperature field can be decomposed into the SHW 
basis functions: 

AT "Jo"-"- 

J-1 

+ ^ ^ ^ 7m,i,i-i/'m,j,l(a;i) (1) 

where Xj^^i and ym,j,i are the approximation and detail 
wavelet coefficients respectively. The first term in (|^) cor- 
responds to a smoother image of the original map, whereas 
the detail coefficients encode the differences between the 
smoothed map and the original. The index j runs over the 
different scales, with J being the resolution of the origi- 
nal map and jo the coarsest resolution considered. Uj is the 
number of pixels at resolution j. Finally, the index m corre- 
sponds to the number of different wavelet functions at each 



Figure 1. In the top figure, the coadded 53-1-90 GHz COBE- 
DMR map in HEALPix pixehsation at resolution J = 7 with the 
customised Galactic cut is plotted. The bottom map shows the 
twelve pixels of the base-resolution of HEALPix (j = 1) and how 
they are affected by such a cut. 

scale required in order to form a complete orthogonal basis 
set. In particular, for a square partitioning, such as that in 
the HEALPix pixelisation, we require three different wavelet 
basis functions at each scale j, and therefore we have three 
different kind of wavelet coefficients 7i,j,i, 72.j.i and 73,j,i. A 
more detailed description of the SHW transform is given in 
Appendix |^. 

As pointed out by MHL, we note that the asymmetry of 
orthogonal wavelet basis functions means that the wavelet 
decomposition is sensitive to the orientation of the input sig- 
nal. Thus any statistics based on the corresponding wavelet 
coefficients are also sensitive to the orientation of the anal- 
ysed signal. 

2.2 Application to COBE-DMR data 

In this paper, we analyse the 4-year COBE-DMR data in 
HEALPix pixelisation. A detailed description of the process 
of map- making is given by Banday et al. (2000). HEALPix 
is an equal area, iso-latitude and hierarchical pixelisation of 
the sphere. The base-resolution comprises twelve pixels in 
three rings around the poles and equator (see Fig. |l|). The 
resolution level of the grid is expressed by the parameter 
A'sidc, that indicates the number of divisions along the side 
of the base-resolution pixel that is needed to reach a de- 
sired high-resolution partition. For convenience we will use 
instead the index j to refer to the scale; this relates to N^ide 
as ATaido = 2-'-\ 

In the case of the COBE-DMR maps, the total number 
of pixels is 49152 (i.e. J = 7, or equivalently A'sidc ~ 64), 
with a pixel linear size of ~ 55'. Therefore, each resolution 
level j corresponds to a scale 55' x 2''~-' , containing nj = 
12 X 4-'"^ pixels of that size. 

We have used in our analysis the coadded 53 A, 53B, 
90A and 90B map (each pixel weighted according to the 
inverse of its noise variance) with the customised Galactic 
cut, which is plotted in Fig. |^ (Banday et al. 2000, Banday 
et al. 1997, Bennett et al. 1996). As mentioned above, we are 
thus using approximately two-thirds of the data as opposed 
to previous analyses that kept only one-third (Face and 
Face 5) of the QuadCube COBE-DMR data. In addition, 
we avoid any possible projection effects in going from the 
sphere to the cube. 

2.3 The non-Gaussianity test 

In our wavelet analysis we have considered separately the de- 
tail coefficients corresponding to m = 1, 2, 3 at each scale j. 
For each value of j and m, we use the corresponding detail 
coefficients "fm,j,i to estimate the skewness S and (excess) 
kurtosis K of the parent distribution from which the coeffi- 
cients were drawn. We therefore obtain the skewness S{j, m) 
and kurtosis K{j, m) spectra for the image. At each value 
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{j,m), the skewness and kurtosis of the parent distribution 
of the wavelet coefficients is given by 

C / 3/2 



K = K4/«:2 



(2) 
(3) 



where k„ is the nth cumulant of the distribution. Follow- 
ing Hobson, Jones & Lasenby (1999) we use fe-statistics (see 
Kenney & Keeping 1954; Stuart & Ord 1994) to obtain un- 
biased estimators of the cumulants, which are then used to 
estimate the skewness S and (excess) kurtosis K. 

In addition to the skewness and kurtosis spectra, we 
measure the scale-scale correlation between the wavelet co- 
efficients at different scales using the following estimator: 



'm,j,{l/i]^m,j + l,l 



m,j,[l/i] L-il = 



(4) 



where [ ] denotes the integer part and rij = 12 x 4-'^^ is the 
number of pixels at the resolution level j. The scale-scale 
correlation coefficient CP{j,m) measures the correlation be- 
tween each type of detail coefficient 7 (i.e. for m = 1,2,3) 
in the two consecutive scales j and j + 1. We restrict our 
analysis to the case p — 2. 

In our non-Gaussianity test, we first obtain the skew- 
ness, kurtosis and scale-scale correlation spectra for the 
coadded 53-1-90 GHz COBE-DMR map. Those wavelet co- 
efficients sensitive to pixels inside the Galactic cut are not 
used in the computations. We then generate 10000 realisa- 
tions of CMB all-sky maps, smoothed with a 7° Gaussian 
beam, and add levels of Gaussian noise corresponding to the 
COBE-DMR data. The CMB maps are drawn from an in- 
flationary/CDM model with parameters i},n — 1, ~ 0, 
h = 0.5, n = 1 and Qrms-ps = 18^K, but we do not expect 
our analysis to be very sensitive to a difference choice of pa- 
rameters (see MHL). We then compute our estimators for 
each of the CDM simulations to obtain approximate prob- 
ability distributions for the S{j,m), K{j,m) and C^{j,m) 
statistics for a CMB signal derived from the chosen model. 
By comparing these distributions with the values obtained 
from the COBE-DMR data, we can obtain the probability 
that our data are derived from a Gaussian distribution char- 
acterised by the chosen power spectrum. 

For each kind of detail coefficient (m = 1, 2, 3), we cal- 
culate the skewness and kurtosis at five different scales, and 
also calculate four scale-scale correlation coefficients. There- 
fore we compute a total of 42 different statistics, which must 
be taken into account when assigning a statistical signifi- 
cance to a given detection. In addition, as we have already 
discussed, this wavelet decomposition is sensitive to the ori- 
entation of the data. Therefore, we perform the entire test 
for three different orientations of the COBE-DMR map. 



3 RESULTS 

The computed S{j,m), K{j,m) and C^{j,m) spectra are 
plotted in Fig. |^ for three different orientations. We rotate 
the 53+90 COBE-DMR map around an axis passing through 
the North and South Galactic poles. Due to the characteris- 
tics of the HEALPix pixelisation, a rotation in this direction 
of 90 degrees simply shifts the twelve base-pixels into each 
other, and thus recovers the original result. Therefore, we 



show the results computed from the COBE-DMR data for a 
rotation of (orientation A, solid squares), 30 (orientation 
B, crosses) and 60 (orientation C, stars) degrees with respect 
to the original orientation of the 53-1-90 COBE-DMR map. 
The open circle corresponds to the average over 10000 CDM 
simulations generated as explained before and the error bars 
indicate the 68, 95 and 99 per cent confidence levels of the 
distributions for orientation A. Note that for convenience 
the K{j, m) and S{j, m) spectra have been normalised at 
each value of (j, m) such that the variance of the distribu- 
tion obtained from the 10000 CDM realisations is equal to 
unity. We do not plot the corresponding error bars for orien- 
tations B and C in Fig. ^ for the sake of clarity and because 
the conclusions derived from the plot remain unchanged. In 
fact, due to the isotropy of the CMB, we expect these error 
distributions to be independent of the chosen orientation, 
although small variations do appear (see Fig. ^ due to the 
presence of anisotropic noise and as a consequence of having 
a pixelised map. 

We see from Fig. ^ that for orientation A (solid squares) 
all the points of skewness and scale-scale correlation spec- 
tra lie within their respective Gaussian probability distribu- 
tions. However, we find a detection of non-Gaussianity at 
a confidence level > 99 per cent for the kurtosis at j — 2, 
m = 3. On the other hand, for orientations B (crosses) and C 
(stars) all the COBE-DMR values are consistent with being 
derived from a parent Gaussian distribution. 

We must be careful when assessing the significance of 
the non-Gaussian detection for orientation A. We must take 
into account that we have computed a total of 3 x 42 = 126 
different statistics and most of them show no evidence of 
non-Gaussianity (see Bromley & Tegmark 1999). Following 
MHL, since the different statistics are not independent, we 
must use Monte-Carlo simulations to estimate the proba- 
bility of having at least one detection of non-Gaussianity 
when computing these 126 statistics at a confidence level 
> 99 per cent even in the case of an underlying Gaussian 
signal. We find that this occurs in 69 per cent of the cases. 
Therefore this analysis does not provide strong evidence of 
non-Gaussianity, in agreement with the results obtained by 
MHL using planar wavelets. 

We may wonder, however, if our conclusions are affected 
by choosing a different set of orientations for the COBE- 
DMR maps. As an illustration of the sensitivity of the com- 
puted statistics to the orientation of the data. Fig. ^ shows 
S{j,m),K{j,m) and C^{j,m) for m = 2 and two different 
scales versus the angle that the COBE-DMR data have been 
rotated around an axis passing through the North and South 
Galactic poles with respect to the original signal. The aver- 
age value and the 68, 95 and 99 per cent confidence levels 
obtained from the corresponding Gaussian distributions are 
also plotted for orientations A,B and C. Note that, in this 
case, the skewness and kurtosis spectra has not been nor- 
malised as in Fig. || It can be seen that the COBE-DMR 
values of the skewness, kurtosis and scale-scale correlation 
at each scale oscillate around the mean obtained from the 
CDM realisations and with a deviation in agreement with 
the confidence limits showed for orientations A,B and C. We 
have also checked that, as expected, we obtain very similar 
variations when rotations are performed around different di- 
rections. Therefore, although the particular numerical values 
of the considered statistics are not rotationally invariant, the 
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Figure 2. The skcwncss, kurtosis and scale-scale correlation spectra for the 53+90 GHz COBE-DMR map are plotted for each type of 
detail coefficient m = 1 (first column), m = 2 (second column) and m = 3 (third column). Solid squares, crosses and stars correspond 
respectively to the the data rotated 0, 30 and 60 degrees around an axis passing through the North and South Galactic poles with respect 
to the input signal (orientations A, B and C). The error bars denote the 68, 95 and 99 per cent limits of the corresponding distribution. 
For convenience, S{j,m) and K(j,m) have been normalised at each value of {j, m) such that the variance of the distribution obtained 
from the 10000 CDM realisations is equal to unity 



general conclusions concerning non-Gaussianity of the data, 
at least in our case, do not seem to depend on the chosen 
orientation. 



4 CONCLUSIONS 

We presented an analysis of the 4- year COBE-DMR data 
(in HEALPix pixelisation) based on spherical Haar wavelets 
(SHW) in order to look for large scale non-Gaussianity of 
the GMB. This analysis is performed using all the available 
data lying outside the Galactic cut. This constitutes about 
two-thirds of the pixels, as compared to just one-third in 



former works that use only Face and Face 5 of the Quad- 
Cube pixelisation. We take into account the sensitivity of the 
method to the orientation of the original signal and present 
the results for three different orientations of the data. We 
also find that the choice of a different set of orientations, 
at least in our case, does not alter the general conclusions 
regarding non-Gaussianity. 

We find that the value of the kurtosis for j = 2, m = 3 
for one of the chosen orientations lies outside the 99 per 
cent confidence level derived from 10000 CDM/inflationary 
realisations, whereas the rest of the statistics show no ev- 
idence of non-Gaussianity. However, since we have a total 
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Figure 3. The S{j, m),K{j, m) and C^{j, m) spectra are shown for values of j = 3, 5 and m = 2 versus the angle that the COBE-DMR 
data have been rotated around an axis passing through the North and South Galactic poles with respect to the original signal. The 
average value and the 68, 95 and 99 per cent confidence levels obtained from the corresponding Gaussian distributions are also plotted 
for values of the orientation angle equal to 0, 30 and 60 degrees (orientations A, B and C, respectively). Note that the skewness and 
kurtosis spectra have not been normalised as in Fig. ti 



of 126 different statistics, the probability of finding at least 
one of them falling outside the 99 per cent confidence level 
even in the case of an underlying Gaussian field is as high 
as 0.69. Therefore, we conclude that an analysis based on 
SHW of the 4- year COBE-DMR data show no evidence for 
non- Gaussianity. 
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APPENDIX A: SPERICAL HAAR WAVELETS 
FOR HEALPix PIXELISATION 

In this Appendix, we give an heuristic approach to decom- 
posing a map in HEALPix pixelisation in terms of spheri- 
cal Haar wavelet coefficients. More mathematical approaches 
are given elsewhere (Sweldens 1995). 

HEALPix is an equal area, iso-latitude and hierarchi- 
cal pixelisation of the sphere. The resolution level of the 
grid is given by the parameter j (or equivalently A^'side, 
A'aide ~ 2-'"^). A level j comprises a total of Uj = 12 x 4-'"^ 
pixels, each of them with equal area fij . Each pixel / at reso- 
lution j, Sj^i, is divided into four pixels Sj+i.ig, Sj+ij, at 
resolution j + 1. In particular, the COBE-DMR resolution 
corresponds to J = 7. 

In order to perform a SHW decomposition of the 
COBE-DMR data, we require one scaling ipjj and three 
wavelet functions ipm.i.i at each scale j and position /. These 
are given by 



1, if xeSj^i 
0, otherwise 



Vj + l.io +^3 + 1,12 ~ (Vj + l.il +^3 + 1,13) 



1 



2/ij+i 
1 



2m. 



'j+i 



(Al) 

(A2) 
(A3) 
(A4) 



where lo,li,l2, h are the four pixels at resolution level j + 1 
contained within the pixel I at level j, and we have taken into 
account that all pixels have equal area at a given resolution. 
The orthogonality of the wavelet functions follows imme- 
diately from their vanishing integral (by construction) and 
their support. The temperature field AT/T on the sphere 



Figure Al. The COBE-DMR map in HEALPix pixelisation at 
resolution J = 7 and the wavelet coefficients at j = 6 are shown. 
Starting from the coefficients at the highest resolution J, the ap- 
proximation and detail coefficients at J — 1 are obt aine d as linear 
combinations of the original image (as shown in Fig A2 ) . The pro- 
cess is repeated over the map formed by the Aj„i j coefficients 
and so on down to the lowest resolution considered jo. 
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Figure A2. The procedure to obtain the wavelet coefficients at 
a resolution level j is shown. Each coefficient is obtained from a 
linear combination of four approximation coefficients at the level 
J + 1. The new Xj,l are just the average over the four original 
ones. 



The value of the g. 



m,j,l, 



factors used to obtain the details is 



— /^j+i or as indicated in the figure. 



can then be written as 



3 J-1 

+ X] X! X! '>m,j,iV'mj,;(2:«) 



(A5) 



where Xj^^i and 7mj,i are the approximation and detail 
wavelet coefficients respectively. The index j runs over the 
different scales with J the resolution of the original map and 
Jo the coarsest resolution considered. 

The proced ure to obtain the wavelet coefficients is as 
follows (see Fig Al). We start with a map at a resolution 
J, identifying each pixel with the Ajj coefhcients. Each ap- 
proximation and detail coefficient at position I and resolu- 
tion level J — 1 is obtained as a linear combination of the 
four corresponding pixels lo,li,l2M at resolution J: 



i=0 

3 

7m,J-l,i = ^ ] gm.J-l,li\j,. 



(A6) 
(A7) 



The approximation coefficients at the new resolution level 
are just the average over the four original coefficients and, 
therefore, correspond to a lower resolution version of the ini- 
tial image. Complementary to those, the details encode the 
difference between the degraded and the original maps. The 
5mj,i coefficients take values ~fij-i or 0, depend- 
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ing on the indices, as shown in Fig. |A^. These coefficients 
are obtained as an integral over the dual function of il)m,j,i 
in the corresponding support (see e.g. Tenorio et al. 1999). 
Therefore, at level J — 1 we have three different sets of de- 
tails and one set of approximation coefficients, each of them 
with nj-i = nj/4 pixels. 

Then, we apply the same procedure to the degraded 
map formed by Aj_i,i , obtaining again a coarser image and 
a set of details corresponding to the resolution level J — 2. 
Repeating this process down to the lowest resolution con- 
sidered in our analysis, jo, we are left with the details for 
resolutions J— 1 to jo as well as the approximation Ajp coeffi- 
cients. With this scheme it is also straightforward to identify 
those coefficients which are descendants of pixels contami- 
nated by the Galaxy and discard them when performing the 
analysis. 

To recover the map we just need to invert the process, 
starting with the approximation and detail coefficients at the 
lowest resolution jo up to the map at the initial resolution 
J. 
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